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SIMULATION  OF  ELECTROSTATIC  MODES  IN  A  MAGNETOPLASMA 
WITH  TRANSVERSE  INHOMOGENEOUS  ELECTRIC  FIELD 

I.  INTRODUCTION 

The  ion-cyclotron  instability  has  been  important  to  both  space  and 
laboratory  plasmas.  In  most  studies,  field-aligned  currents  or  ion  beams 

have  served  as  the  driving  mechanism. ^  Recently  some  laboratory 

2  i 

experiments"’"'  have  reported  ion  cyclotron  instability  in  circumstances 

where  neither  field-aligned  currents  nor  ion  beams  provide  a  satisfactory 

free  energy  source,  and  the  existence  of  an  inhomogeneous  electric  field  is 

conjectured  to  play  a  role.  Furthermore,  unstable  ion  cyclotron  waves  have 

4 

been  reported  in  connection  with  double  layers  and  shocks  in  spa<_e 
plasmas."*  In  all  these  cases  a  localized  electric  field  perpendicular  to 
the  ambient  magnetic  field  is  an  intrinsic  feature  of  the  plasma 
equilibrium. 

In  order  to  understand  the  role  of  the  d.c.  electric  field  in  the 
excitation  of  the  ion  waves  around  the  ion  cyclotron  frequency,  Ganguli  et 
al.^  developed  a  nonlocal  kinetic  theory.  In  the  early  version  of  the 
theory  the  gradients  of  the  electric  field  were  ignored  so  as  to  avoid  the 
confusion  with  the  Kelvin-Helmholtz  (K-H)  instability  which  is  explicitly 
dependent  on  the  second  derivative  of  the  electric  field.  They  suggested  a 
possible  mechanism  based  on  the  coupling  of  positive  and  negative  energy 
waves,  in  a  plasma  equilibrium  containing  a  uniform  magnetic  field  along 
with  a  transverse  localized  electric  field,  that  can  excite  the  ion  waves. 
Dependent  on  the  parameters  such  as  the  magnitude  of  the  d.c.  electric 
field,  these  waves  can  have  real  frequencies  around  the  ion  cyclotron 
frequency. 

Recently  the  kinetic  theory  has  been  generalized  by  Ganguli  et  al^  to 
include  continuous  gradients  in  the  d.c.  electric  field,  and  in  the  weak 

t  t 

shear  limit  (i.e.  Vg(x)/S3.  <<  1,  where  is  the  spatial  derivative  of  the 
E  X  B  drift  and  Q  is  the  ion  cyclotron  frequency)  the  distinctions  between 
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the  new  ion  waves  and  the  kinetic  K-H  waves  have  been  studied.  They 

obtained  a  general  dispersion  relation  at  the  second  order  differential 

equation  level  (see  the  Eq.  (29)  in  Ganguli  et  al.^).  For  the  localized 

2 

transverse  electric  field  profile  (~  sech  x)  and  the  initial  density  as 

shown  in  Fig.  1,  the  growth  rates  of  the  kinetic  K-H  mode  (Fig.  2  (a))  were 

obtained  for  various  values  of  u  =  k^k^  where  k(|  and  k^  are  the  components 

of  the  wave  vector  in  the  z  (parallel  to  the  magnetic  field  Bq)  and  the  y 

2  2 

directions.  The  K-H  mode  is  localized  below  b  =  k  p.  =  1  and  k  L  -  1, 

y  i  y 

where  L  is  the  scale  length  associated  with  the  d.c.  electric  field  and  p. 
=  v^/S^  is  the  ion  gyro-radius.  Furthermore,  the  growth  rates  are  reduced 
strongly  as  u  increases.  The  calculated  growth  rates  of  the  ion-cyclotron¬ 
like  modes  with  u  =  0.038  are  plotted  in  Fig.  2(b).  This  mode  becomes 

unstable  where  b  >>  1  and  k  L  >>  1.  The  unstable  domains  for  both  modes 

y 

are  clearly  distinctive. 

In  this  paper  we  investigate  the  new  kinetic  ion-cyclotron-like  ( I — C— 
like)  instability^  driven  by  a  localized  electric  field  perpendicular  to 
the  ambient  magnetic  field  by  means  of  particle  simulations.  The  empnasis 
in  this  paper  is  to  demonstrate  the  existence  of  the  new  mode^  and 
distinguish  it  from  the  well  known  Kelvin-Helmholtz  branch  and  the  kinetic 
modification  of  this  branch.  The  initial  loading  of  particles  in  phase 
space  is  another  significant  issue  in  regard  to  the  simulation  of  this 
instability,  because  the  most  obvious  choice  of  initial  distribution 
function,  i.e.,  a  Maxwellian  shifted  in  the  v  velocity  component  due  to  the 

g 

E  x  B  drift,  is  not  a  self  consistent  equilibrium.  Previous  attempts  at 
simulating  the  lC-like  modes  have  encountered  difficulties  as  a  result  of 
this  fact.  Therefore,  we  will  also  attempt  to  clarify  this  issue.  In  Sec 
II  the  simulation  model  and  result'  are  presented  and  concluding  remarks 


and  discussions  are  described  in  Sec.  III. 
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II.  SIMULATION  MODEL  AND  RESULTS 


A  two-dimensional  electrostatic  code  is  used  which  retains  the  full 
dynamics  of  the  ions  in  three  dimensional  velocity  space.  Electrons  are 
treated  by  the  guiding  center  approximation  in  the  perpendicular  direction 
while  the  parallel  motion  is  treated  exactly.  We  use  a  system  length  L^  = 
128A,  Ly  =  32 A  or  64A,  where  A  is  the  grid  size  which  is  equal  to  the 

_  2 

electron  Debye  length,  X  and  n  X  =  36.  The  external  electric  field  is 

0  0  0 

2 

applied  in  the  form  of  Eq;<(x)  =  E^sech  [  ( x  -  64)/L]  in  the  x-direction 

which  produces  E  X  B  drift  in  the  y-direction  given  by  V_(x)  =  -E  (x)/B  . 

E  ox  o 

For  the  present  problem,  in  which  a  nonuniform  d.c.  electric  field 
transverse  to  a  uniform  magnetic  field  is  present,  the  initial  distribution 
function  constructed  out  of  the  constants  of  motion  was  provided  by  Ganguli 
et  al . , 7 


foa,H)  =  N  exp(-(3H)g(0 


(1) 


where  N  =  no(3/2rt) 


2 

1/v  ,v  is  the  thermal  velocity  and 


g(0  =  exp  {gleTai/m  +  VE2(  l)/2  J }  n(  0"1/Z 


where  U(  O  =  1+Vg(  £,)/$?,  [|  =  d[j/d£,,  £  =  x  +  (v  £,))/$  is  the  guiding 

2  2  2 

center  position,  H  =  (v^  +v^  iv^  )/2  *  eTfxl'ni  is  the  total  eneigy  and 

T(x)  is  the  electrostatic  potential  associated  with  the  d.c.  electric 
field.  For  the  purpose  of  initial  loading  in  a  computer  simulation,  the 


3 


distribution  (1)  can  be  expressed  in  terms  of  the  physical  position  x  (for 
details  see  Ganguli  et  al.2), 


2nf 


3exp{-|[v2  + 


ol 


(vy-VE(X))' 

_ rM) _ 


n 


[1  +  0(e)), 


/noo 


(3) 


It  is  interesting  to  note  that  the  distortion  in  the  gyro-orbit  introduced 
by  the  sheared  transverse  d.c.  electric  field  leads  to  a  sustainable 
difference  in  the  temperature  in  the  two  dimensions  transverse  to  the 

t 

uniform  magnetic  field.  For  weak  shears  (i.e.  V£  /Q  <<  1),  (3)  reduces  to 


2nfol  =  P{1  ♦ 


ev  E(x> 
To 


l(vy-VE(x))2  -  v2)}  exp(-  §(v2+(v  -VE(x))2))  •  (4) 


We  see  that  for  weak  shear  the  use  of  a  shifted  Maxwellian  for  the 
initial  loading  is  acceptable  and  the  system  will  adjust  to  the  above 

t 

distribution  automatically  with  only  minor  relaxation.  If  the  V„  /0 
correction  is  ignored,  then  (4)  becomes  identical  to  a  Maxwellian  shifted 
by  the  magnitude  of  the  E  X  B  drift  in  the  y  component  of  the  velocity. 

O 

Such  an  initial  distribution  was  used  by  Pritchett  and  Coroniti  and 

9  ' 

Pritchett  but  they  extended  it  to  higher  shear  values  (V^  /Q  >  1  so  that  r> 

<  0)  where  this  initial  distribution  is  no  longer  acceptable.  If  the 

initial  state  is  far  from  equilibrium  then  there  will  be  additional 

artificial  free  energy  which  may  lead  to  strong  relaxation  of  the  initial 

nonequilibrium  state.  Such  strong  relaxation  from  a  nonequilibrium  initial 

condition  is  invariably  accompanied  by  a  substantial  release  of  free 

energy,  which  in  turn  leads  to  a  noisy  dynamic  state  quite  different  from 
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the  quiescent  initial  state  essential  for  simulation  of  the  instability. 

/ 

Therefore,  for  larger  values  of  Vg/Q,  the  interpretations  of  the  simulation 

8  9  8 

results  in  these  papers  ’  (especially  the  conclusion  of  Pritchett  ,  that 

the  Kelvin-Helmholtz  instability  will  always  dominate  over  the  new 

instability^)  remain  dubious. 

In  the  following  we  shall  ignore  the  order  e  terms  but  retain  the  shear 

corrections  and  employ  (3)  for  the  initial  loading  of  our  simulations. 

Later  on  we  shall  also  use  the  shifted  Maxwellian  scheme  for  a  weak  and  a 

moderately  strong  shear  case  and  compare  the  results  with  the  more 

appropriate  loading  scheme  given  in  (3). 

We  first  choose  B  /B  (=k,./k  )  =  0.0075  and  L  =  64X.  .  The  minimum  b 
oy  o  ll  y  y  e 

becomes  0.21.  As  shown  in  Fig.  2(a),  the  kinetic  K-H  mode  has  maximum 
growth  rate  around  b  =  0.2.  Thus  we  can  expect  the  simultaneous  excitation 
of  both  the  branches. 

Other  parameters  used  are  9  /oo  =4,  T./T  =  3.5,  m./m  =  100,  L  = 

p  e  pe  ’  l  e  l  e  ’ 

11A,  a  =  E"  /4itn  T  =  0.2,  V°  =  0.59v.  where  n  ,  V°,  and  v.  are  the 
Oxee  E  it  e  E  it 

averaged  electron  density,  the  peak  value  of  Vg(x)  and  the  ion  thermal 
velocity,  respectively.  The  time  evolutions  and  their  spectra  are  shown 
for  the  modes  (0,  3)  which  corresponds  to  b  1.89,  k^p.  =  1.37  and  k^L  = 

3.23,  and  (0,4)  which  corresponds  to  b  =  3.38,  k^p^.  =  1.84  and  k^L  =  4.31, 
in  Fig.  3.  As  shown  in  Fig.  3. (a)  the  mode  (0,  3)  begins  to  grow  clearly 
and  emerge  from  the  background  thermal  noise  around  9. t  =  150,  and  goes 
into  the  nonlinear  stage  around  9.  t  =  300.  The  detailed  wave  analysis 
shows  that  in  the  linear  stage  this  mode  has  real  frequency  around  0.5ST. 
As  shown  in  Fig.  3  ( b ) ,  the  peak  of  the  frequency  spectrum  is  located 

around  0.29..  The  frequency  of  this  mode  becomes  lower  in  the  nonlinear 
stage.  This  is  a  common  feature  for  all  the  modes.  The  real  frequency  of 
the  (0,  4)  mode  saturates  around  0.59.  as  shown  in  Fig.  3  (d).  Also  from 
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fig.  3  (b)  and  (d)  we  find  a  number  of  smaller  peaks  around  higher 
harmonics. 

An  estimate  of  the  growth  rates  of  several  modes  in  order  to  identify 
the  fastest  growing  mode  indicates  that  the  mode  (0,  4)  which  corresponds 
to  b  =  3.38,  has  the  maximum  growth  rate  y/'B.  -  0.025.  In  Fig.  4  we  plot 
the  estimated  growth  rates  of  the  (0,1)  to  (0,5)  modes  against  b.  It 
should  be  remarked  that  the  mode  (0,  1)  which  corresponds  to  b  =  0.21  and 
kyL  =  1,  and  falls  on  the  K-H  branch  (the  solid  curve  in  Fig.  4)  has  a 
smaller  growth  rate  (y/Q.  ~  0.015)  than  the  higher  modes  which  form  the  new 
branch.  Clearly  it  is  the  ion-cyclotron-like  branch  with  the  mode  (0,4) 
that  actually  dominates.  The  maximum  growth  rate  of  this  ion-cyclotron- 
like  mode  is  located  around  large  b  (=  3.38  and  kyL~4)  and  large  real 

frequency  (~  ,  in  the  linear  stage).  Thus,  unambiguous  distinction  is 

made  between  the  kinetic  K-H  mode  and  the  new  ion-cyclotron-like  mode. 
Besides,  we  have  identified  both  the  branches  of  the  system  in  our 
simulation  and  for  the  given  set  of  parameters  we  find  that  the  dominant 
mode  does  not  belong  to  the  K-H  branch;  instead  it  lies  on  the  new  I-C-like 
branch  which  has  higher  frequency  and  shorter  wavelength. 

The  average  ion  flow  velocities  v^(x)  are  shown  in  Fig.  5.  In  the 
earlier  linear  stage  (B.t  <  80)  there  is  hardly  any  reduction  in  the  peak 
of  the  velocity  profile  or  any  relaxation  in  the  topology  of  the  profile. 
In  the  linear  stage  we  find  that  the  real  frequency  is  around  the  ion- 
cyclotron  frequency  and  its  higher  harmonics.  This  is  in  sharp  contrast  to 

g 

the  simulations  of  Prichett  ,  where  there  was  a  strong  relaxation  of  the 
initial  profile  right  from  the  very  onset  of  the  simulation.  The  reason 
for  this  strong  relaxation  can  be  traced  to  the  employment  of  the 

O 

simplified  initial  loading  scheme  of  a  shifted  Maxwellian  in  Prichett  ; 
which  is  not  acceptable  when  the  magnitude  of  Vg/Q.  is  significant.  Ue 
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discuss  the  critical  importance  of  the  initial  loading  scheme  later  on. 

In  the  late  quasilinear  (S.t  =  250)  and  the  nonlinear  (S^t  >  300) 
stages  the  wave  amplitudes  become  large  and  can  interact  with  and  reduce 
the  magnitude  of  the  d.c.  electric  field.  Thus  in  the  final  stages  of  the 
simulation  we  see  a  reduction  in  the  magnitude  of  Vp(x).  Also  since  the 

real  frequency  is  proportional  to  k  V  (x)  (see  Ganguli  et  al.^),  we  find 

y  t, 

that  in  the  quasilinear  stage  the  real  frequency  starts  reducing  along  with 
the  magnitude  of  V„(x).  Further  reduction  in  the  real  frequency  occurs  in 
the  fully  nonlinear  phase. 

In  the  previous  simulation  we  included  a  slight  density  gradient.  The 

profile  of  transverse  electric  field  E  (x)/E  ,  the  initial  density 

ox  ox’  J 

profiles  for  the  electrons  and  the  ions  are  shown  in  fig.  6(a),  (b),  and 
(c).  The  density  ratio  of  the  maximum  to  minimum  n^/n^,  equals  to  1.1. 
The  density  gradients  are  localized  at  the  half  value  of  the  transverse 
electric  field  as  shown  in  Fig.  6.  The  reason  that  the  ion  density  profile 
is  not  as  smooth  as  the  electron  density  profile  is  due  to  the  fact  that 
the  ions  have  larger  gyro-radius.  The  location  of  the  guiding  centers  of 
ions  are  exactly  the  same  as  those  of  the  electrons.  However,  the  gyro- 
motion  of  ions  smears  out  the  density  gradients.  This  particular  profile 
will  be  used  for  all  the  simulations  in  this  article  (except  Fig.  7(a)  and 
(b)). 

In  the  preceeding  discussions  we  have  established  the  existence  of  the 
I-C-like  branch  and  distinguished  it  from  the  K-H  branch.  Ue  have  also 
demonstrated  that  at  least  for  some  parametets  the  I-C-like  branch  can 

dominate.  Hence,  in  order  to  economize  on  the  CPU  time  we  shorten  L  to 

v 

32A  in  the  remainder  of  this  paper  to  isolate  the  I-C-like  branch  and  study 
the  characteristic  features  of  this  branch  only  (except  Figs.  9  and  11). 


Note  that  the  mode  (0,4)  will  now  correspond  to  the  mode  (0,2)  and  so  on. 
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The  time  evolutions  of  the  perturbations  of  the  electrostatic  potential 
for  the  mode  (0,2)  with  and  without  density  gradients  and  their  spectra  are 
shown  in  Fig.  7.  In  these  cases  we  use  B^/B^  =  0.005,  keeping  the  rest  of 
the  parameters  unchanged.  The  case  with  the  density  gradients  has  the 
growing  wave  with  slightly  larger  growth  rate,  larger  amplitude,  and  lower 
real  frequency,  as  shown  in  Fig.  7.  Inclusion  of  even  a  slight  density 
gradient  shows  enhancement  in  the  growth  rate  and  the  amplitude.  This  is 
in  contrast  to  the  K-H  instability  where  density  gradient  leads  to 
reduction  of  the  growth  rates  of  the  instability^.  The  fact  that  density 
gradients  lead  to  stronger  growth  for  the  I-C-like  modes  was  predicted  by 
the  theory^. 

The  dependence  of  the  real  frequency  on  Vg(x)  is  now  examined.  The 

time  evolutions  of  the  perturbations  of  the  electrostatic  potential  for  the 

mode  (0,2)  corresponding  to  V°  =  0.59v.  and  2.0v.  and  their  spectra  are 

h  it  it 

shown  in  Fig.  8.  Here  B  /B  =  0.005  and  the  rest  of  the  parameters  remain 

oy  o  r 

unchanged.  the  real  frequencies  in  the  linear  stages  (not  necessarily  the 
largest  peak)  are  found  around  0.3Q.  (b)  and  1.0C.(d).  As  expected  from 

the  theory^,  we  find  that  the  real  frequencies  of  the  instability  are 
approximately  proportional  to  kyV°. 

The  diffusion  of  the  guiding  center  of  ions  in  the  case  of  =  64A  is 
shown  at  R..  t  =  0,  80,  160,  240,  and  320  in  Fi0.  9.  A  subset  of  ions  with 
initial  position  such  that  58  <  x  <  64,  are  followed  in  time,  and  the  local 
guiding  centers  of  ions  are  plotted  at  each  time.  At  the  time  ( R.  t:  =  160) 
when  the  waves  begin  to  grow  out  of  the  thermal  noise,  slight  waving  of  the 
column  is  observed.  At  the  beginning  of  the  nonlinear  stage  (R.t=240),  the 
column  is  twisted.  At  the  nonlinear  stage  the  ions  are  strongly  diffused 
due  to  the  vortices. 
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The  evolution  of  vortices  at  the  nonlinear  stage  (Q. t  -  240)  is  shown 
in  the  case  of  =  32A  in  Fig.  10.  The  electrostatic  potential  in  the  x-y 
plane  has  strong  vortices  at  the  center.  The  strong  vortex  moves  along  the 
direction  of  the  E  x  B  drift.  As  shown  in  Fig.  10,  the  excited  waves  are 
well  localized  around  the  center.  In  the  case  of  L  =  64A  similar  strong 

y 

vortex  is  also  created  which  moves  along  the  direction  of  the  E  x  B  drift. 

In  the  later  nonlinear  stage  the  number  of  the  vortices  is  increased. 

The  local  electric  fields  E  at  x  =  69,  y  =  14  in  the  case  of  L  = 

x,y  y 

64A  is  diagnosed,  and  the  power  spectrum  is  obtained.  As  shown  in  Fig. 
11(a),  around  t  =  300  strong  sharp  peaks  are  produced  by  the  passing 

vortices.  The  spectrum  has  several  peaks  around  Gk  which  correspond  to  the 
real  frequencies  in  the  linear  stage.  the  largest  peak  is  produced  by  the 
vortices  in  the  nonlinear  stage. 

The  profiles  of  overall  wave  energy  in  the  case  of  Ly  =  32A  are  shown 

2 

in  fig.  12.  The  time-averaged  excited  wave  energy  (E^)  is  localized  around 

the  center  of  the  system.  Obviously,  the  excited  wave  energy  is  much 

larger  than  the  thermal  energy  outside  the  electric  field.  The  time- 

2 

averaged  polarization  electric  field  energy  with  minus  sign  -(Ex)  '  is  also 
plotted,  which  is  much  smaller  than  the  excited  wave  energy. 

Now  we  demonstrate  the  importance  of  the  proper  initial  loading.  We 
first  examine  the  weak  shear  case  where  V°/LS2  =  0.24,  and  compare  the 

difference  in  results  produced  by  the  simplified  and  the  improved  initial 
conditions.  As  expected  and  noted  earlier  there  is  no  remarkable 
difference  in  the  weak  shear  limit.  However,  as  shown  in  Fig.  13  (a)  and 
(c)  the  growth  rate  of  the  mode  (0,  2)  with  the  improved  initial  conditions 
is  slightly  larger  and  the  amplitudes  the  nonlinear  stage  are  slightly 


reduced.  Furthermore,  at  r  =  8o  the  a-eiage  ion  I  low  velocities  witli  our 
improved  initial  condition  remain  almost  unchanged  whereas  the  use  of  the 
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shifted  Maxwellian  initial  condition  leads  to  relaxation  in  the  velocity 
profile,  although  marginal  in  this  case  since  the  shear  strength  is  weak. 
Apparently,  our  initial  condition  improves  the  simulation  results  even  in 
the  weak  shear  limit. 

The  simulations  with  the  larger  transverse  electric  field  (V^  =  2-Ov.^ 
which  corresponds  to  a  =  2.24)  with  density  gradients  have  been  carried  out 
keeping  the  other  parameters  the  same  as  in  Fig.  13.  Now  the  shear  of  the 
d.c.  electric  field  is  moderate  (V°/LS.  =  0.86).  Therefore,  we  expect 

Ci  1 

greater  improvements  in  the  simulation  results  by  employing  our  initial 
conditions  provided  in  Eq.  (3)  over  the  simple  initial  condition  (shifted 
Maxwellian).  Perturbations  of  electrostatic  potential  for  the  mode  (0,2) 
corresponding  to  the  two  different  initial  conditions  are  compared.  As 
shown  in  Fig.  14,  the  mode  with  our  improved  initial  conditions  has  much 
larger  growth  rate  and  larger  amplitude  at  the  end  of  the  simulation.  This 
demonstrates  that  with  the  improved  initial  condition  not  only  does  the 
initial  noise  level  decreases  but  also  the  instability  becomes  much  more 
coherent  and  pronounced. 

The  evolutions  of  the  average  ion  flow  velocities  v^(x)  corresponding 
to  the  two  different  initial  conditions  discussed  earlier  are  also  compared 
for  the  larger  transverse  electric  field  case.  As  expected,  with  our 
initial  conditions  the  early  relaxation  and  reduction  of  the  average  ion 
flow  until  SJ. t  =  80,  was  greatly  reduced  as  (see  Fig.  15(b)).  In  this 
case,  the  free  energy  is  solely  dedicated  to  excite  the  large  amplitude 
waves  as  shown  in  Fig.  14(c).  The  amplitudes  of  the  excited  waves  become 
quite  large  at  around  S2. t  =  240  and  therefore  can  interact  with  the  d.c. 
electric  field  to  reduce  it's  magnitude.  At  the  same  time  the  excited 
strong  vortices  cause  large  diffusions  of  ions  (see  Fig.  9(e))  and 
consequently  the  average  ion  flow  velocities  are  greatly  reduced. 
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III.  DISCUSSIONS 


We  have  investigated  the  new  ion-cyclotron-like  waves  driven  by  a 
localized  transverse  electric  field  by  means  of  a  numerical  simulation  with 
the  assistance  of  the  nonlocal  kinetic  theory.  The  linear  theory  shows 
that  the  growth  rates  of  the  kinetic  Kel vin-Helmhol t z  modes  are  strongly 
reduced  with  increasing  u,  and  they  become  unstable  only  where  b  <  1  and 
kyL  as  1.  On  the  other  hand,  the  new  modes  have  larger  real  frequencies  and 
become  unstable  at  larger  b  >  1  and  k^L  >>  1. 

We  have  performed  a  number  of  simulations  for  the  new  kinetic  ion 

modes.  Results  show  that  ion-cyclotron-like  waves  are  excited  in  the 

regions  where  the  E  X  B  drift  is  localized.  The  linear  growth  rates  of 

several  modes  are  estimated  from  the  wave  analysis  of  the  simulation.  This 

linear  analysis  shows  that  the  (0,  4)  mode  which  corresponds  to  large  b  and 

large  k^L,  has  the  maximum  growth  rate.  Clearly,  these  are  not  K-H  modes 

(See  Fig.  (?)).  Further,  the  simulation  results  show  that  density 

gradients  help  to  enhance  the  growth  rates.  This  in  in  contrast  to  the 

established  properties  of  the  Kelvin-Helmhol tz  instability^.  However, 

like  the  K  H  mode  the  real  frequencies  of  this  instability  are 

approximately  proportional  to  k  V°.  This  is  predicted  by  the  nonlocal 

y  h 

theory . 

Furthermore,  we  established  the  critical  importance  of  the  initial 
condition  which  includes  the  modification  due  to  the  gradient  of  the  d.c. 
electric  field.  ^  While  this  improved  initial  condition  is  preferable  even 
for  weak  shears  of  the  d.c.  electric  field,  it  becomes  crucial  when  the 
shear  is  large.  We  have  demonstrated  that  for  moderate  shear  values,  out 
initial  condition  provides  t lie  quiescent  initial  state  essential  for 
simulations  of  the  instability.  TIkm -of  m  c: ,  r  hr-  initial  relaxation  of  t  lie 
average  ion  flow  velocity,  which  is  unavoidable  when  the  simple  initial 
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condition  (the  shifted  Maxwellian)  is  used  for  shear  values  that  are 
significant,  is  greatly  reduced. 

The  nonlinear  phenomena  such  as  diffusion  and  coalescence  of  vortices 
are  investigated.  In  the  linear  stage,  smaller  vortices  are  generated  and 
larger  vortices  with  the  lower  real  frequencies  are  dominant  in  the 
nonlinear  stage.  In  the  nonlinear  stage  the  ions  diffuse  strongly  due  to 
the  large  scale  vortices. 
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The  time  evolutions  ((a),  (c)  and  their  power  spectra  (b),  (d)>  of 
the  perturbations  of  electrostatic  potential  for  two  modes  in  the 
case  of  Ly  =  64A  with  the  slight  density  gradients. 

(a),  (b):  k|(/Xe(|  =  0.0075  x  2n  x  3/64  =  2.21  x  103,  kyAe||  =  2 it  x 

3/64  =  0.295. 
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theoretically  predicated  Kelvin-Helmholtz  branch  for  identical 


The  time  evolution  of  perturbations  of  electrostatic  potential  of 
the  (0,  2)  mode  in  the  case  of  L  =  32A  ((a),  ( c ) )  and  their  pover 
spectra  (b),  (d)>,  (V°  =  0.59  V.^  Boy/Bo  =  0.005). 

(a),  (b):  no  density  gradient 
(c),  (d):  with  the  gradients. 
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Fig.  13.  Comparison  of  the  results  with  two  different  initial  conditions 

with  density  gradients  as  shown  in  Fig.  5  (The  mode  (0,2)  in  the 

case  of  L  =  326, B  /B  =  0.005). 
y  py  o 

(a),(b):  the  simplified  initial  condition. 

(c),(d):  the  improved  initial  condition. 
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Fig.  14.  The  time  evolution  of  perturbation  of  electrostatic  potential 

((a),  (c))  and  their  pover  spectrum  ((b),  (d))  in  the  case  of  a  = 

E2  Man  T  =  2.24,  v°  =  2.0v.„,kllX  ..  =  0.005  x  2n  x  2/32  =  1.96  x 
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(a),(b):  the  simplified  initial  conditions. 

(c),(d):  the  improved  initial  conditions. 
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(a):  the  simplified  initial  conditions. 
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